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Meta-Orbital Transition in Heavy-Fermion Systems: Analysis by Dynamical Mean 
Field Theory and Self-Consistent Renormalization Theory of Orbital Fluctuations 
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We investigate a two-orbital Anderson lattice model with Ising orbital intersite exchange 
interactions on the basis of a dynamical mean field theory combined with the static mean 
field approximation of intersite orbital interactions. Focusing on Ce-based heavy-fermion com- 
pounds, we examine the orbital crossover between two orbital states, when the total f-electron 
number per site m is ~ 1. We show that a "meta-orbital" transition, at which the occupancy of 
two orbitals changes steeply, occurs when the hybridization between the ground-state f-electron 
orbital and conduction electrons is smaller than that between the excited f-electron orbital and 
conduction electrons at low pressures. Near the meta-orbital critical end point, orbital fluc- 
tuations are enhanced and couple with charge fluctuations. A critical theory of meta-orbital 
fluctuations is also developed by applying the self-consistent renormalization theory of itinerant 
electron magnetism to orbital fluctuations. The critical end point, first-order transition, and 
crossover are described within Gaussian approximations of orbital fluctuations. We discuss the 
relevance of our results to CeAi2, CeCu2Si2, CeCu2Ge2, and related compounds, which all have 
low-lying crystalline-electric-field excited states. 

KEYWORDS: dynamical mean field theory, two-orbital Anderson lattice model, self-consistent renormal- 
ization theory, critical end point, orbital fluctuations 



1. Introduction 

After the discovery of the unconventional supercon- 
ductivity in CeCu2Si2, 1 various heavy-fermion com- 
pounds have attracted great attention. They exhibit 
many interesting phenomena, 2 such as quantum criti- 
cal behaviors associated with the quantum critical point 
(QCP) of magnetic phase transitions under pressure, 
magnetic field and chemical substitutions, non-Fermi liq- 
uid properties near the QCP, and various types of un- 
conventional superconductivity. It is expected that the 
Cooper pairs in such heavy-fermion superconductors are 
mediated by magnetic fluctuations, since the supercon- 
ductivity always occurs near the magnetic phase. 

Over the past ten years, the possible nonmagnetic- 
fluctuation-mediated superconductivity at high pressures 
separated from the superconductivity near the mag- 
netic QCP has attracted much attention in some heavy- 
fermion compounds, such as CeCu2Si2, CeCu2Ge2, and 
related materials, 3-7 since magnetic fluctuations are not 
developed there. 8 Near the pressure where the supercon- 
ducting transition temperature has a peak as a function 
of pressure, electric resistivity shows a linear tempera- 
ture T dependence, with residual resistivity being en- 
hanced. 5-7 As pressure increases, Kadowaki- Woods ratio 
gradually varies from the value for strongly correlated 
systems to that for weakly correlated systems. 9 

To explain superconductivity without magnetic fluc- 
tuations, Onishi and Miyake 10 proposed that Ce-valence 
fluctuations lead to unconventional superconductivity at 
high pressures in these systems on the basis of an ex- 
tended periodic Anderson model (Ex-PAM). The ex- 



istence of the valence transition in the Ex-PAM was 
confirmed by more elaborate numerical calculations 
later. 11 ' 12 The properties near the critical end point 
(CEP) of valence fluctuations explain the resistivity 7, 13 
and changes in Kadowaki- Woods ratio. 9 

Experimentally, when the f-electron charge per site 
varies, owing to the electrostatic potential change around 
a Ce site, the lattice constant should also vary corre- 
spondingly, and some of the phonon branches should be 
affected by this change. However, no direct evidence for 
the steep changes in Ce valence has been reported. 

In this paper, we will focus on the crystalline-electric- 
field (CEF), i.e., orbital states in these systems. The ma- 
terials mentioned above all have low-energy excited CEF 
states. CEF excitations are observed by inelastic neutron 
scattering experiments. 14, 15 Resistivity data also indi- 
cate the existence of excited CEF states, showing a two- 
peak structure as a function of temperature T. The two 
peaks correspond to two different Kondo temperatures. 
More importantly, the two peaks merge as pressure in- 
creases, and the pressure where the two peaks merge ap- 
proximately coincides with the pressure where supercon- 
ducting transition temperature increases as mentioned 
above. 5 The two-peak structure and its variation as a 
function of pressure are well explained by the analysis of 
an orbital-degenerate Anderson lattice model. 16 Interest- 
ingly, these behaviors are also observed in CeAl2, which 
is a prototype of heavy-fermion compounds, at around 3 
GPa, while superconductivity has not been observed. 17 
A neutron scattering experiment revealed that there are 
also low-energy CEF states in CeAk). 18 
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Thus, it is natural to consider that the orbital (CEF) 
fluctuations in these systems play an important role in 
realizing the unconventional superconductivity at high 
pressures. In this paper, we will examine the orbital vari- 
ations as a function of pressure on the basis of a dynami- 
cal mean field theory (DMFT) 19 and show that a "meta- 
orbital" transition or crossover occurs when hybridiza- 
tions between conduction and f-electrons depend on the 
orbital and are sufficiently small compared with the f- 
clcctron energy level. We also find that the orbital fluc- 
tuations couple with f-electron charge degrees of freedom, 
leading to the variation in f-elcctron occupancy when the 
orbital crossover or transition occurs. We will also try to 
construct an effective orbital fluctuation theory. 

This paper is organized as follows. In §2, a two-orbital 
Anderson lattice model is analyzed on the basis of a dy- 
namical mean field theory with Wilson's numerical renor- 
malization group method used as the impurity solver. A 
theory of critical fluctuations near the critical end point 
of orbital fluctuations is developed in §3 on the basis 
of the self-consistent renormalization theory. Finally, §4 
shows the summary of the present paper. 

2. Two-orbital Anderson Lattice Model with 
Ising Anisotropic Orbital Interaction 

In this section, we investigate orbital fluctuations in a 
two-orbital Anderson lattice model with Ising orbital in- 
tersite exchange interactions on the basis of the DMFT 19 
combined with the static mean field approximation of in- 
tersite exchange interactions. As a solver for the effective 
impurity problem in the DMFT, we use Wilson's numeri- 
cal renormalization group (NRG) method, which is pow- 
erful for investigating zero-temperature properties of the 
impurity problem. We will show the zero-temperature 
properties of this system, including the variations in the 
orbital occupancy, the f-electron density of states, and 
the zero-temperature phase diagram. 

2.1 Model 

We investigate the periodic Anderson lattice model on 
a Bethe lattice with two localized f-orbitals (a = A or B) 
hybridizing with a single conduction electron band. 16 We 
also include Ising orbital-orbital intersite interactions in 
our model and the Hamiltonian is given as 

H = ~ ^n)4* c j° + [McLf™ + h - c -) 

ija iaa 

+ {e f 
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Here, Ci a and f iaa represent the annihilation operators 
for conduction and f-electrons at site i, spin a, and 
orbital a, respectively. The f-electron number and the 
orbital operator are defined aS Tlifau = fiaa fiaa and 
T iz = 12a( n ifAo- - n if Ba), respectively. J* is the Ising 



orbital interaction between the nearest neighbor sites, 
which is assumed as > 0, and z n . n . is the number 
of nearest-neighbor sites. Since, in the DMFT, intersite 
correlations are neglected, we introduce this orbital inter- 
action phenomenologically in order to take into account 
the correlations. Orbital interactions are generated by 
a fourth-order perturbation with respect to v a , starting 
from the localized limit in a two-orbital Anderson lat- 
tice model, as shown in Appendix A. Hamiltonian (1) 
includes only one of the many interactions generated in 
fourth-order perturbations, such as spin-spin exchange 
and interactions between spin-orbital degrees. Since we 
are interested in the orbital fluctuations in a paramag- 
netic phase, we restrict ourselves to take into account 
the Ising orbital interaction in this paper. The notations 
used for other parameters are conventional. 

Throughout this section, we set the nearest-neighbor 
hopping t = 1/2, which corresponds to half of the band- 
width D = 2t = 1, and we concentrate on analyzing 
zero-temperature T = properties of Hamiltonian (1) 
and fix the total electron number per site n = 1.8, and 
SfA — —0.5 in the unit of D. As for the value of J*, we 
will mainly discuss the case of = 0.02 while briefly 
showing the results for J„ = as a reference. Since we 
are interested in Ce-based heavy-fermion compounds, we 
analyze the f 1 configuration, i.e., the configuration with 



rif 



y 



1. This is the reason why spin ex- 



change interactions between A- and B- orbitals are ab- 
sent in Hamiltonian (1), since they are irrelevant in the 
f 1 configuration. In the following, U a and U' are set to 
be infinite in order to project out f^-f 4 configurations. 

2.2 Analysis by the dynamical mean field theory 

In the DMFT, a lattice model is mapped to an effec- 
tive impurity one with an effective conduction electron 
bath. 19 In our case, the effective impurity model is an im- 
purity Anderson model with three localized orbitals hy- 
bridizing with single-band conduction electrons. In this 
paper, the intersite orbital interaction is approximated 
in the mean field level with a uniform amplitude. The 
effective impurity model is given as 

HcS = [ha\ a a k a + (^La\ a c a + h.c.) - fij^ 4 c <x 

ka * a 

+ [ Va (4/«<T + h - C + foe. - tiflafav] 



^ U a n fat n fa i + U' ^2 nfAanfBa' 
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Here, Sf a includes the contributions of the mean field 
value of T tz : e fA = e fA -J^(T z ) ande /B = e fB + J^(T z ), 
where (• • • ) indicates the expectation value and (T z ) — 
{Ti z ). The operators at the impurity site are defined as 
fo a a = faa and co<r = c a . The effective bath is repre- 
sented by the bath conduction electron a\ a hybridizing 
only with c a . The dispersion and hybridization of the 
bath conduction electrons are characterized by and 
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\>k, respectively, which are calculated self-consistently in 
the DMFT. Although one can integrate both aha and Co- 
in a path-integral formulation, it is useful to employ the 
form (2) for NRG calculations. 

To solve the impurity problem in the DMFT, we em- 
ploy Wilson's NRG method, 20,21 which is powerful for 
investigating the low-energy properties of the impurity 
problem. To calculate one-particle Green's functions cor- 
rectly, we use the full-density-matrix NRG developed re- 
cently 22,23 In the Bethe lattice, the self-consistent equa- 
tion of the DMFT reads 



A(w) = — V = -G c 

h 



(3) 



where G c {w) = -i dte- luJt (T{ CrT (t), 4(0)}) being the 



' dte- wt {T{c a {t),cU 
local retarded Green's function for conduction electrons. 
In the NRG calculations in the DMFT self-consistent 
loops, we use the logarithmic descretization parameter 
A = 1.8 and mainly keep 300 low-energy states at each 
NRG iteration. 20 We also confirm that the results do not 
change when we increase the number of retained states 
up to 700. 

2.2.1 Meta-orbital transition 

In this subsection, we will show the numerical results 
of DMFT+NRG calculations. We will demonstrate that 
a meta-orbital transition and crossover occur in this sys- 
tem. We will also demonstrate that the orbital degrees 
of freedom couple with the f-electron charge ones. 

Figure 1 shows the EfB dependence of orbital occu- 
pancy (T z ) for seven sets of vb with fixed va — 2ug. 
As EfB decreases, (T z ) decreases from ~ 1 to ~ —1, in- 
dicating that the dominant orbital occupation changes 
from the A-orbital to the B-orbital. The variation in 
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Fig. 1. (Color online) (T z ) vs EfB for several values of vg = 
0.511,4, £fA = —0.5, and J* = 0.0 and 0.02, and the filling is 
fixed to be n = 1.8. The arrow indicates the typical variation 
in (T 2 ) for = 0.02 as a function of pressure for Ce-based 
heavy- fermion systems. 
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Fig. 2. (Color online) (a) rif vs EfB and (b) n c vs EfB for several 
values of vb = 0.5^. The other parameters are the same as 
those in Fig. 1. 



(T z ) as a function of EfB is gradual when Vb{= 0.5?;^) is 
large, while it becomes steeper as vb decreases. Finite 
also enhances this steep change. Orbital occupancies for 
Jo = 0.0 with vb =0.1 (squares) and vb — 0.09 (filled 
squares) change smoothly, while those for = 0.02 with 
the same parameters show much steeper changes. 

A CEP is located between vb — 0.095 and 0.1 for 
Jo = 0.02. For vb = 0.09, a first-order transition occurs. 
One can see a clear hysteresis, as shown in Fig. 1. For 
vb = 0.095, although we have not detected a hysteresis 
in our calculation, it seems that there is a first-order 
transition. 

These behaviors are analogues of meta-magnetism. In 
the present system, e/b and (T z ) correspond to the mag- 
netic field and magnetization in the language of meta- 
magnetism, respectively. This meta-orbital transition is 
expected to occur in systems with several low-energy 
CEF states as the pressure increases. When the param- 
eters are tuned, the CEP is realized. In Fig. 1, we also 
show the expected variations in vb, £fB, and (T z ) in 
Ce-based compounds, with increasing pressure, indicated 
by an arrow for = 0.02. Note that, as the pressure 
increases, both the hybridization and f-clectron energy 
levels are expected to increase in Ce-based heavy-fermion 
compounds. 

Figure 2 shows the e/b dependences of the total f- 
electron number rif and the conduction electron number 
n c with fixed n = rif + n c = 1.8. Note that, for vb < 0.1, 
nf is close to 1 for a small EfB, while nf ~ 0.9-0.95 for a 
large EfB- This is because vb < va- When EfB is suffi- 
ciently small, for example, EfB — —0.55 and vb = 0.095, 
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(T z ) ~ —0.9, as shown in Fig. 1, which indicates that the 
f-electron occupies essentially the B-orbital. On the other 
hand, when ejb is large, the f-electron occupies mainly 
the j4-orbital. In these two limiting cases, it is sufficient 
to consider a one-orbital Anderson lattice model. We ob- 
tain a low Kondo temperature Tk that leads to rif ~ 1 for 
the small EfB with the relevant hybridization vg < va, 
while a high Tr that leads to rif < 1 for the large EfB 
with the relevant hybridization va > vb- 

As for the symmetry, the orbital, f-electron number, 
and conduction electron number belong to the same ir- 
reducible representation (scalar). Thus, as shown in Fig. 
2, when the meta-orbital transition or crossover occurs, 
the f-electron number changes accordingly. We note that, 
across the transition and the crossover, the f-electron 
number is still close to 1, when va and vb are not very 
large as expected in heavy- fermion systems. This point 
contrasts the results of the valence transition driven by 
the Coulomb repulsion between the f- and conduction 
electrons, Uf c , in the Ex-PAM, where a large variation 
in f-electron number occurs, owing to the large Uf c as- 
sumed. 10-12 

2.2.2 Variation in density of states 

Figure 3 shows the density of states (DOS) of f- 
electrons —\m.G f a,b(w) / it for vb = 0.1 = 0.5va and 
several values of EfB- Here, Gjaab)^) represents the 
local A(B)-orbital f-electron retarded Green's function 
with the energy uj. In addition to the lower Hubbard- like 
peak located at the energy uj ~ £fA,B ~ —0.5, there are 
several peaks in the low-energy region, reflecting the CEF 
excited states. The low-energy peak near uj ~ 0.2 in the 
DOS of the B-orbital for s/b — —0.44 moves to a lower 
energy as SfB decreases, while a peak appears at a posi- 
tive energy in the DOS of the A-orbital as EfB decreases, 
reflecting the meta-orbital crossover discussed previously. 
There is a sharp peak near uj = in — Ivs\G f a{b) I ^ ■ 
One can also see that there is a hybridization gap like 
structure for the A(£?)-orbital DOS when ejb is larger 
(smaller) than approximately EfB ~ —0.53, as shown in 
Figs. 3(b) and 3(c). 

To investigate the low-energy DOS in more detail, we 
show the DOS at uj — as a function of EfB in Fig. 4, 
where we define Gtot = G fA + G fB + G c . Since ImG c 
does not show notable changes, we do not show it. The 
DOS of the ^4-orbital at uj = is large when EfB is large, 
while it is small when s fs is small. Correspondingly, the 
DOS of the £?-orbital at uj = is large for a small EfB 
and small for a large EfB- In Fig. 4, one can also see a 
clear hysteresis in the DOS at uj = for vb = 0.09. 

If we employ the variations in hybridization and f- 
electron level as pressure increases as shown by the arrow 
in Fig. 1, we observe that the DOS at the Fermi level de- 
creases as pressure increases. The steepness of this vari- 
ation depends on the distance to the CEP from the pa- 
rameter path that corresponds to the physical pressure. 




Fig. 3. (Color online) (a) f-electron density of states vs w for 
vg = 0.1 and va = 0.2 for several values of £/b- Full and 
dashed lines represent the DOS for the A- and B-orbitals, re- 
spectively. The values of ejg decrease from bottom to top as 
e fB = -0.44, -0.45, -0.46, -0.47, -0.48, -0.49, -0.50, -0.505, 
-0.51, -0.515, -0.52, -0.525, -0.528, -0.53, -0.532, -0.535, 
-0.538, -0.542, -0.547, -0.555, and -0.56. Each line is shifted 
by 0.3. The other parameters are the same as those in Fig. 1. 
(b) Low-energy structure of A-orbital density of states, (c) Low- 
energy structure of B-orbital density of states. 

2.2.3 Phase diagram 

To summarize this section, we show the schematic 
ground state phase diagram of the two-orbital Ander- 
son model with two different hybridizations (vb < va) 
and ferro-orbital intersite interactions in Fig. 5. As dis- 
cussed above, we observe a first-order transition line and 
its CEP. In Fig. 5, the "A(B)-rich" region indicates the 
region with a larger f-electron occupancy in the A(B)- 
orbital. The "A- rich" region is located in the larger EfB 
region, where the effective mass is smaller than that in 
the "B-rich" region, owing to the large va > vb- We also 
call it "light" Fermi liquid instead of heavy Fermi liquid. 
This light Fermi liquid can be regarded as the state with 
larger valence fluctuations than the B-rich state, since 
va > vb- Near the CEP, f-elcctron charge fluctuations 
and also conduction electron ones are enhanced. The ef- 
fective mass also changes steeply near the CEP, reflecting 
the variation in orbital occupancy. 

Since our analysis is based on the DMFT and the 
mean field approximation for intersite orbital interac- 
tions J*, we overestimate the effects of J%. As shown 
in Fig. 1, the EfB dependence of (T z ) for J I — is much 
weaker than that for J* = 0.02. We have not carried out 
DMFT+NRG calculations for a very small v B < 0.09, 
since very small hybridizations make the NRG calcula- 
tion difficult. When — > 0, whether the first-order tran- 
sition occurs or not is beyond the scope of the present 
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Fig. 4. (Color online) f-electron density of states at the Fermi 
level (w = 0) vs £fg. ( a ) Total density of states — ImC?tot(0)/7r, 
(b) —lmGfj\(0)/ir, and (c) — ImG f g (0) / it . The other parameters 
are the same as those in Fig. 1. 



analysis. However, it is noted that when one of the hy- 
bridizations is zero, there must be a first-order transition 
even when J* = 0. To clarify whether this first-order 
transition exists in the region where the hybridization 
is finite, more elaborate numerical calculations that can 
take into account the intersite correlations are needed. 

As for magnetic instabilities, it is expected that mag- 
netic phases will appear in the phase diagram, especially 
in the region of the smaller vb and £/s, corresponding 
to the low-pressure region in typical Ce-based heavy- 
fcrmion systems. Determining the phase diagram includ- 
ing the magnetic phases and their magnetic structures is 
beyond the scope of the present paper and we leave it 
as a future study. It is also interesting to investigate the 
case that the magnetic transition and the orbital CEP 
occur simultaneously. Such a situation might lead to a 
new universality class of magnetic phase transitions. 

3. Self-Consistent Renormalization Theory of 
Nonmagnetic Orbital Fluctuations 

In this section, we will develop a self-consistent renor- 
malization theory for the orbital fluctuations discussed in 
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Fig. 5. (Color online) Schematic ground state phase diagram in 
vg-Efg plane for vg = aVj\ with fixed constant a < 1. Arrows 
indicate examples of the parameter changes as the physical pres- 
sure P increases. The dotted line indicates the orbital crossover 
line where (T z ) ~ 0. The critical end point is represented by a 
filled circle and the first-order transition line is indicted by a solid 
line. The hysteresis related to the first-order transition is not in- 
dicated. The left side of the first-order transition corresponds 
to the state with a larger effective mass due to the smaller hy- 
bridization vg , while the right side of the first-order transition 
corresponds to that with a smaller effective mass, denoted as 
light Fermi liquid. 



§2 on the basis of the DMFT. Since the DMFT does not 
take into account the long wavelength fluctuations, we 
use the one-loop self-consistent approximation scheme of 
the fluctuations known as the self-consistent renormal- 
ization (SCR) theory developed by Moriya and Kawa- 
bata for itinerant electron magnetism. 24-26 When we ap- 
ply the SCR theory to the nonmagnetic orbital fluc- 
tuations, we need to take into account the tempera- 
ture dependence of the order parameter itself and the 
temperature dependence originating from the order pa- 
rameter in the mode- coupling corrections. Our discus- 
sion in this section leads to results similar to those of 
the SCR theory for meta-magnetic fluctuations 27 and 
Kondo- volume-collapse transition in heavy-fermion sys- 
tems. 28 It is noted that our formulation naturally de- 
scribes both the order parameter and the mode-coupling 
corrections in an equal footing based on a variational 
principle. We will discuss the orbital fluctuations in 
three-dimensional (d = 3) systems in this paper, al- 
though the application to the cases with d ^ 3 is straight- 
forward. 

3.1 Landau- Gintzburg- Wilson action 

The Landau-Gintzburg- Wilson action for nonmagnetic 
fluctuations generally has odd order terms with respect 
to field variables, which is given as 



S 



hcj) 
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Here, <p p is the nonmagnetic field with (f>* = <p- p , and 
p represents the momentum and energy. N and T are 
the number of lattice sites and temperature, respec- 
tively. In terms of S, the free energy F is given by 
cxp(— F/T) = J T>(j>exp(—S). We neglect the momen- 
tum and energy dependence of the coefficients b and c, 
which are irrelevant from simple power counting in the 
context of the renormalization group. For the stability of 
the system, c > 0. h is the field conjugate to the uniform 
nonmagnetic field <po; thus, does not depend on p. Since 
we are interested in the orbital transition without break- 
ing translational symmetry, we use the Ornstein-Zernike 
form of a p as in the ferromagnetic fluctuation: 24-26 



a p = a + A\q\ 2 + C|w„|/|q| 



where q and u> n are the momentum and the Bosonic Mat- 
subara frequency, respectively. Existence of the odd order 
terms in eq. (4) leads to the first-order transition, CEP, 
and also crossover, as shown in Fig. 5. 

3.2 Modified SCR theory 

In the SCR theory, the non-Gaussian terms in the 
free energy F are approximated by mean field decou- 
pling. 24 ~ 26 - 29 This approximation corresponds to varia- 
tionally determining the coefficients of the Gaussian fluc- 
tuations. In our case with the presence of h and b terms, 
we need to take into account the order parameter itself 
in addition to mode-coupling terms. We will show that 
a modified Gaussian action can describe the first-order 
transition, the CEP, and the crossover. 

To describe the order parameter variations within the 
Gaussian approximation, we assume the following varia- 
tional action: 

SeS = \ X! r P$P$-p + + 7^o</>o, (6) 



p#0 

Op + 5, 



(7) 



where 6 and h are variational parameters. On the basis 
of Gibbs-Bogoliubov-Feynman inequality, they are deter- 
mined by minimizing f2: 30 



Q = O e flf + T(S — S e ff), 

fieff = -Tl0gZ e S, 



^eff 



(S) 



X>0exp(-S* cff ), 



Vcj)S exp(-S eS )/Z eS . 



(8) 
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(10) 
(11) 
(12) 



Defining the order parameter (f> as 4> = y/T/Nfo, we 



obtain 

n 



, = + (h-h)(cb)- S -(cb 2 )+b(^)+c^) 
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Since the effective action S e s is Gaussian, we can analyt- 
ically calculate {<p n ) with n being any integers, and they 
are given as 

_h_ 
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Differentiating fl c g with respect to S, we obtain 
h-h + 36-^ — + 3b (0 2 ) - 8c(0) 3 - 6b(<jy} 2 



dS 



+ 



d_ 

ds 



p#0 



1 



Using eq. (15) and defining 

p ' 

we obtain 

h-h + 3bX- 3b((j)) 2 - 8c(0) J 



(18) 



(19) 



+ 



--+6cX + 3b(<j))+6c(<j)) 2 



d_ 

85 



+ X 



= 0. 
(20) 

Note that when h = h = b = {<p) = 0, cq. (20) leads to 
the self-consistent equation 

5 = 12cX, (21) 

which coincides with the conventional expression. 24,26 
Similarly, differentiating n e g with respect to h leads to 

d(<f>) 



h-h + ZbX -3b((py -8c(0) c 



6cX 



dh 

36(0) +6c<^>> 2 l4-(0 2 > =0.(22) 
J oh 
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Combining eqs. (20) and (22), we obtain 

h = h + 3bX - 3b((f)) 2 - 8c((p) 3 , (23) 

5 = \2cX + Qb{(p) + 12c(0) 2 . (24) 

Equations (23) and (24) are the self-consistent equations 
in our theory. 

Now, we derive the equation of state for (</>). Using eqs. 
(14), (23), and (24), the equation of state is derived as 

h + 3bX + (ao + 12cX){4>) + 3b(<P) 2 + 4c(</>) 3 = 0. (25) 



(a) (a - a^.)l2T A = -0.0001 



(b) o = fl 0c 



Parameterizing (<p) as {<p) = . 
poses, eq. (25) is rewritten as 



bj (4c) for later pur- 



apb 
\~c 



2b 3 



(4c) 2 



36^ 

4c 



Equation (26) is the basic equation used to determine 
the order parameter S(f>. Note that Sip appears only in 
the form of Sip 2 in X from the definition of Sep and eq. 
(24), and the effect of fluctuations on h in eq. (25) is 
canceled out by the shift of — > Sip — 6/ (4c). The fluc- 
tuation term appears only in the coefficient of Sep in eq. 
(26). This equation describes the first-order transition, 
the CEP, and also the crossover while taking into ac- 
count the mode-coupling correction X within the Gaus- 
sian fluctuation theory. A similar equation of state was 
introduced in ref. 28, where the renormalization of h was 
neglected in eq. (25). 

Figure 6 shows the temperature and h dependences 
of the order parameter Sep. Here, the critical value arj c 
is given by eq. (36) and the parameters To, T A , and x c 
introduced in Appendix B are fixed as To = 20 K, Ta = 
50 K and x c = 300. In Fig. 6(a), there is a first-order 
transition for h/2T A ~ 0.05749 below T/T - 0.01. In 
Fig. 6(b), the parameter ao is fixed as ao = ao c ; and thus, 
at zero temperature, the CEP is realized at h ~ 0.05476. 
In Figs. 6(c) and 6(d), there is no phase transition and 
Sep varies smoothly as functions of T and h. 

3.3 Zero-temperature critical end point 

In this subsection, we investigate the nonmagnetic or- 
bital fluctuations near the zero-temperature CEP in de- 
tail. Details of the calculations are summarized in Ap- 
pendix B. 

Near the zero-temperature CEP, ro ~ 0, X in eq. (19) 
is expanded in terms of tq as 



X = K +K iro +JC(T), 



(27) 



with 



ro 



ao 



\2c\Ko + K x ro 



3b 2 



YlcS 



1 



1 - YlcKx 



3b 2 

a - — + UcScp 2 + l2c(K Q + JC(T)) 
4c 



I II 



hflT A 
(c) (Of,- a c)/2Z> 



h/2T A 
(d)Ob-a 0c )/27i = 



! 



+ 12cX)5(P + 4c5(p 3 = 0. (26) 



h/2T A 



h/2T 4 



Fig. 6. (Color online) Contour plot of order parameter 8<f> in 
temperature-ft plane for b/2T A = -0.01, c/2T A = 0.01, T = 20 
K, T A = 50 K, and x c = 300. (a) (o - a 0c )/2T A = -0.0001, (b) 
ao = a c, (c) (a -a 0c )/2T A = 0.0001, and (d) (a -a 0c )/2T A = 
0.001. In (a), a finite temperature critical end point appears at 
h/2T A ~ 0.05749 and T/T ~ 0.01, below which a first-order 
transition line exists. In (b), a zero- temperature critical end point 
exists at h/2T A ~ 0.05746, where the order parameter <50 = 
for T > 0. In (c) and (d), changes in the order parameter are 
smooth, indicating the crossover. 



As is well known, /C(T) oc T 4 / 3 in three-dimensional sys- 
tems with the dynamical exponent z — 3 24 > 26 (see Ap- 
pendix B). Now, substituting eqs. (27) and (28) into eq. 
(26), we obtain 



Ho 



Ao + A X K{T) Sep + 4C S(p 3 = 0, 



where 



Ho 



= h- — 



Ao = a — 



Ax 



apb 
4c " 

3^ 
4c 
12c 



26 3 
(4c) 



\2cK-i 



4C 



1 - Y2cK x 

36cKx 



1 - YlcKx 
>0, 



a 



36^ 
~4c 



4c 



1 



(30) 



(31) 
(33) 



(34) 



1 - YlcKx- 

In the case that the Sep dependence of K , Kx, and K.(T) 
is negligible, eq. (30) can be derived by minimizing the 
effective free energy T e ff given as 

1 r 



F e g = H 



A +AxlC(T) 



5<p z + C S(p 4 



(35) 



and 



Ko > 0, Kx < 0, and K,(T) > 0. 



(29) 



which is nothing but the Landau's free energy for ferro- 
magnetism in the presence of the magnetic field —Ho 27 
Since JC(T) vanishes at the zero-temperature CEP 
as T 4 / 3 , the zero-temperature CEP is realized when 
Ho = Ao = 0, which is equivalent to ro = and h = 0. 
( 28 As a result, the critical value of the order parameter is 
Sip = 0, leading to (<p) c = — 6/ (4c). In terms of the orig- 
inal parameters, we obtain the conditions for the zero- 
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temperature CEP as 



a-oc 



K. = 



3/r 2 

4c 
b 3 

- 3bK , 



(4c)= 



(36) 
(37) 



at zero temperature. Using eqs. (36) and (37), eqs. (31) 
and (32) are rewritten as 



Hq — h — h c , 
Aq = (a - a Qc 



YlcKy 



(38) 
(39) 



1 - YlcKx 

Let us discuss the temperature dependence of 5(f> near 
the zero-temperature CEP in the case of Co > as ex- 
pected for a small c > 0. 

First, we consider the case that all the parameters, 
h, do, 6, c, A, and C are constant, i.e., temperature- 
independent. Then, H = and A = in eq. (30), 
when the CEP is realized at zero temperature. We obtain 
5<j> = or 

4i(0) 



8<j>* 



-MT), 



(40) 



4C (0) 

where j4i(0) and Co(0) are the zero-temperature values 
of A x and C , respectively. Since A 1 (0) > 0, /C(T) > 0, 
and Co(0) > 0, the solutions of eq. (40) are pure imagi- 
naries. Thus, 6(f) = for all temperatures without taking 
into account the temperature dependence in the param- 
eters. This corresponds to the order parameter variation 
in Fig. 6(b) for h/2T A ~ 0.05746. 

Secondly, we consider that the leading temperature 
dependence of the parameters is proportional to T 2 as 
expected in Fermi liquid states. In this case, eq. (30) be- 
comes 



Ai{Q)K{T)8<l> + 4C (0)^ 3 ~ -H„T 2 , (41) 

where H ~ HqT 2 . Thus, we obtain 6(j) oc T 2 / 3 , since this 
temperature dependence is self-consistent with JC(T) oc 
T 4 / 3 . Setting K{T) = kT 4 / 3 with n > 0, and 5(f) = ipT 2 / 3 , 
we obtain 



A^iup + 4C (0)<p 3 



(42) 



When Ai(0) > and C o (0) > 0, eq. (42) always has 
only one solution ip, which is either positive or negative, 
depending on the sign of H . The first-order transition 
occurs as H" changes for Ai(0) < and Co(0) > 0. 

It is noted that <f> couples any quantities and thus, in 
general, this nonanalytic temperature dependence 8<p oc 
T 2 / 3 appears in any physical observables near the zero- 
temperature CEP. For example, the magnetic suscepti- 
bility Xs has, rather than the usual T 2 dependence, the 
nonanalytic temperature dependence as 

Xs - 1 (T)^ x ; 1 (0) + A^ + .-- . (43) 

Here, A is the coupling constant between magnetic and 
nonmagnetic fluctuations such as A</>^ p M p M_ p in the 
action, where M p is the magnetic fluctuation at the mo- 



mentum and the energy p. Thus, there might be a finite 
temperature region where Xs(T) oc T~ 2 / 3 is observed, 
when xj 1 (0) is suppressed by some factors. Since the 
present theory is phenomenological, the microscopic dis- 
cussions on this issue are beyond the scope of this paper. 

As for the other thermodynamic and transport quan- 
tities, we obtain the same critical behaviors as those in 
the ferromagnetic QCP with z = 3 and d = 3. The spe- 
cific heat coefficient diverges as C/T <~ — log T, 31 and 
the resistivity is given as p <~ T 5 / 3 , 32 in the asymptotic 
limit. It is noted that, at a higher temperature, p oc T as 
is the case for all the bosonic fluctuations. 

3-4 Critical end point at finite temperature 

Here, we discuss the temperature dependence of 5<p 
near a finite-temperature CEP at T = T* > 0. 

As in the case of the zero-temperature CEP discussed 
in §3.3, we first consider the case that all the parameters 
have no temperature dependence. In this case, 8<j> = for 
T > T* and 5<j) oc y/KjT*) - JC(T) for T < T*. Thus, 
we obtain 5<f> oc (T* - T) 1 / 2 near T ~ T*, exhibiting 
the mean field behaviors. This is understood by noting 
that our theory is essentially mean field approximation 
and that the Ho = line in the parameter space corre- 
sponds to the case of the Ising model without magnetic 
fields, identifying A Q + A^T) as T — T c , where T c is 
the magnetic transition temperature. 

In the second case, the temperature dependence of 
all the parameters is proportional to T 2 as derived in 
Fermi liquid states. Expanding the temperature depen- 
dence from the critical end point temperature T*, we 
obtain 



~ 



dH 



dT 



(T-T*)+4C (T*)5(/) 3 . (44) 



T=T* 



Here, we assume T* <IC Tfl, where Tfl is the so-called 
Fermi liquid temperature below which the Fermi liquid 
behaviors appear. In eq. (44), the temperature depen- 
dence in the second term in eq. (26) or (30) is on the 
higher order in (T—T*)/T*, and thus we neglect it. Then, 
the asymptotic temperature dependence of the order pa- 
rameter S(f> near the finite-temperature CEP is 



S<f> = ± 



1 



4C (T*) dT 



T=T* 



IT* - T\3. 



(45) 



Comparing this with the result at the zero-temperature 
CEP as discussed in Sec. 3.3, one can see that there is 
a classical-quantum crossover from 8<j) oc \T — J 1 *] 1 / 3 to 
5<j) oc T 2 / 3 as T* decreases in the second case. 

Physically, it is unrealistic to consider H = for all 
temperatures in real materials; thus, the second case is 
more realistic for experimental situations. However, the 
first case is pedagogically important, since it clearly ex- 
hibits this transition belonging to the Ising-type univer- 
sality class as expected. It is also important to note that 
8(f) oc T 2 / 3 is realized only when Tfl is well defined, i.e., 
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when the Fermi liquid theory is valid and the condition 
T/T FL <C 1 is realized. 

4. Discussion and Summary 

We have discussed orbital fluctuations in a two-orbital 
Anderson lattice model on the basis of the DMFT in §2 
and constructed a self-consistent critical theory of the 
orbital fluctuations in §3. 

In the first part of this paper in §2, we have investi- 
gated meta-orbital transition in the two-orbital Anderson 
model with anisotropic intersite orbital-orbital interac- 
tions, and the schematic phase diagram is shown in Fig. 
5 as the summary of §2. 

First, our results demonstrate that a meta-orbital 
transition or crossover with steep changes in orbital oc- 
cupancy occurs as pressure increases, only when the hy- 
bridization between the conduction electrons and the f- 
electrons at ambient pressure is sufficiently small. For 
simplicity, let us consider that the pressure affects the 
hybridizations only. Evidently, if the ambient pressure is 
located on the right side of the first-order transition line 
in Fig. 5, the meta-orbital transition and crossover do not 
occur as pressure increases, since vb increases as pressure 
increases. On the other hand, if the ambient pressure is 
located on the left side of the first-order transition line, 
it is possible to realize them, and if the parameters are 
tuned, the CEP is realized. Here, we note that the f- 
electron ground state multiplet on the right side of the 
first-order transition line is the A-orbital, while it is the 
.B-orbital on the left side of the line. Thus, it is con- 
cluded that the meta-orbital transition or crossover is 
realized only when the hybridization between the con- 
duction electrons and the ground state multiplet of f- 
electrons at ambient pressure is smaller than that be- 
tween the conduction electrons and the excited CEF mul- 
tiplet. These observations are valid when the dominant 
effect of the pressure is the variation in va.b- 

As far as we know, there is no compound that shows a 
first-order meta-orbital transition as pressure increases. 
This might indicate that the ferro-orbital interaction is 
small in the real materials or transverse orbital, and other 
interactions smear out the first-order transition. As men- 
tioned in §2.2.3, it is important to clarify whether the 
first-order transition occurs in the model without J*. 
Since J* has been phenomenologically introduced to in- 
corporate the orbital-orbital interactions in the DMFT 
calculations, we need to investigate the J I — model on 
the basis of more elaborate methods that can take into 
account both the local and intersite interactions. 

Another important result shown in §2 is that the meta- 
orbital transition or crossover affects the f-electron occu- 
pancy. As shown in Fig. 2, f-electron occupancy slightly 
changes, when the meta-orbital transition or crossover 
occurs. This is understood by noting the hybridizations 
va 7^ vb- This contrasts with the large valence change 
discussed in the context of the valence transition driven 
by the large Coulomb repulsion Uf c , where the f-electron 



occupation number nf cannot be close to 1, when the 
CEP is realized. The large Uf c also affects the conduc- 
tion electrons. Our recent study of the Ex-PAM based 
on DMFT+NRG calculations reveals that, when the va- 
lence crossover occurs as pressure increases for the large 
Uf c comparable to the bandwidth, the resulting density 
of states of the conduction electrons near the crossover 
point shows a Hubbard-like peak deep below the Fermi 
level and a sharp peak near the Fermi level, which in- 
dicates the existence of strongly correlated conduction 
electrons. 

In the second part of this paper in §3, we have devel- 
oped a self-consistent rcnormalization theory of the or- 
bital fluctuations. There, we have paid special attention 
to the order parameter dependence of the mode-coupling 
term and derived self-consistent equations, which lead 
to an equation of state for the order parameter. We 
have succeeded in describing the first-order transition 
and crossover including the CEP of the orbital fluctua- 
tions. The temperature dependence of the order parame- 
ter oc T 2 / 3 is consistent with that in the Kondo-volumc- 
collapse theory 28 and possibly consistent with that in the 
valence transition, 33 since the symmetries of these order 
parameters are the same. 

As for the temperature dependence of the resistivity, 
it is shown 7 that the valence fluctuations lead to the 
T-linear dependence near the CEP, assuming a small A 
in eq. (5) from the results of the slave boson theory. 10 
It is noted that the T-linear temperature dependence is 
realized in all the bosonic fluctuations at high temper- 
atures. For example, phonon scatterings lead to such a 
behavior at high temperatures. In the case of the orbital 
fluctuations, it is not clear whether the A coefficient in 
eq. (5) is sufficiently small in the present analysis. The 
precise determination of the A coefficient needs a more 
quantitative analysis. When A <~ in eq. (5), we obtain 
JC{T) oc T 2 / 3 near the CEP 33 and the equation of state 
(30) leads to 5(f) cx T 4 / 3 when all the parameters vary 
as a function of T 2 at low temperatures. Note that this 
temperature dependence should be cut off at the scale 
Ta , below which the asymptotic behaviors characteristic 
to d — 3 and z — 3 discussed above appear. 

In our theory, it is expected that CeCu2Si2 and re- 
lated compounds and CeAl2 at ambient pressure are on 
the left side of the line of the first-order transition in Fig. 
5, i.e., in the i?-rich phase. As for superconductivity, it 
is expected that compounds that show superconductiv- 
ity away from their magnetic QCP pass near the CEP 
in Fig. 5 as pressure increases, while nonsuperconduct- 
ing materials, such as CeAl 2 , arc not sufficiently close 
to the CEP, assuming that magnetic phases exist on the 
left side in the phase diagram in Fig. 5. Since there are 
many aspects that affect the realization of superconduc- 
tivity, such as the energy scale of the fluctuations and 
band structures, 34 more quantitative model and analysis 
are needed in order to clarify the details of the supercon- 
ductivity realized in these systems. 
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In summary, we have investigated orbital fluctuations 
in heavy- fermion systems. We have shown that a meta- 
orbital transition occurs as pressure increases and the or- 
bital variations affect the density of states at the Fermi 
level and also f-electron occupancy. We have also con- 
structed an effective critical theory for the orbital fluc- 
tuations, which describes both the first-order transition 
and crossover on the basis of a variational principle. Our 
results give a possible scenario for the superconductivity 
and non-Fermi liquid behaviors associated with orbital 
fluctuations in heavy- fermion systems. 
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Appendix A: Localized limit and intersite inter- 
actions 

In this appendix, we will derive the intersite orbital 
interaction via a fourth-order perturbation in v a from 
the localized limit of the Hamiltonian (I) without the J* 
term. Note that J§ is not exactly the same as the orbital 
intersite interaction that will be derived in this appendix. 
Jl has been introduced phenomcnologically in order to 
take into account the intersite orbital correlation in the 
DMFT calculation in §2. 

Let us consider Hamiltonian (I) without J*. In the 
limit of U a , U' — > oo, and Sf a — > — oo keeping 
U a /\£fp\ 1 and VaVp/sf^ finite, the local excited con- 
figuration of f° contributes to the low-energy properties 
of this system only as virtual states and configurations 
with P(n > 2) are neglected. Thus, the effective Hamil- 
tonian can be written by the localized spin and orbital 
degrees of freedom interacting with conduction electrons 
via exchange interactions: 

Hex = y ~ ^ij) c ia c jcr + ^ ] ^a c Lr( S )cru' c ia l 

ijcr ictcrcr' 

+ A J2 T iz + J2 J ABc\ a {s c )v a ,c ia ,T ix • S 4 

i icrcr' 

+ Kab y ' n icTi Z + K AB y ' ni C Ti X , (A 

i i 

where Si represents the f-electron spin operator at site i 
and the conduction electron spin is given as = cr^/2, 
where <7 M with /i = x, y or z is the Pauli matrix. n ic is the 
conduction electron charge operator n^ c = J^a c \o c i°- As 
for the orbital degrees of freedom of f-electrons, T ix = a x 
and T iz — a z acting on the orbital space at site i and the 
two eigenvalues of T iz are 1 for the A-orbital and —I for 
the B-orbital. The coupling constants are given as 

v 2 

Ja = (A-2) 

- £ f 



Jab 



K 



Va - V% 



- £ f 



ab — 



-4e f 



(A3) 



(A4) 



Here, we assume EfA ~ £/b = £/ and ignore the dif- 
ference between sja and e/s in the denominator of J a , 
Jab, K ab , and K AB . The difference is taken into ac- 
count in A = £fA — £fB- 

Now, intersite interactions in the second-order pertur- 
bations of J a , Jab , K AB , and K AB are obtained straight- 
forwardly. Among them, the Ising-type orbital-orbital in- 
teraction is given as 



(A5) 



*,3 



where \ij is the conduction electron charge susceptibility. 
Note that, when vb ~> va (or va S> Vb), the coupling 
constant Jb (or J a) and \K AB \ is much larger than Jab 
and K AB . In this limit, the spin-spin and Ising orbital 
interactions on the fourth order in v a are more important 
than those in the other sectors, such as the transverse 
orbital and spin-orbital coupled interactions. 

Appendix B: Detailed calculations in SCR the- 
ory 

In this appendix, we summarize the detailed calcu- 
lations of the mode coupling correction X in the SCR 
theory in order to make this paper self-contained. 

Let us parameterize X in eq. (19) as 



X 



T 



1 



N ^ r o + A|q| 2 + CK||q|- ' 



(Bl) 



with 6=1. Following the discussions by Misawa et 
eq. (B-l) is reduced to 



29 



X = JC(0)+JC(T), 



AC(0) = — 
S, n 



dq / duj 
o Jo 



Cujq- 



,1+0-1 



IK j f qc f 
JC(T) = =m / dq duj 

K JO Jo 



[q e (r +Aq 2 )} 2 + (Cujf 
n B (u)Cujq d+0 - 1 



(B-2) 
(B-3) 

[q e {r {) + Aq 2 )] 2 + (Cuj) 2 ' (B ' 4) 

where is the Bose distribution function ns(w) = 
Al{e u,T - 1) and if d =27r d / 2 /[(27r) d r(d/2)], with T being 
the gamma function. In terms of dimensionless parame- 
ters, eqs. (B-3) and (B-4) are rewritten as 

T a d f°°, , , [ Xc a ^ +e - 1 dx 
Ta J 



JC(0) 



[x d (y + x 2 )} 2 + z> 2 ' 



(B-5) 



AC(T) 



2T Q d 
~Ta~ 



zdz 



„d+e-i 



lio [x%y + x 2 )/t] 2 + z 2 



(B-6) 



Here, y = r /2T A = (a + 12cA + 66(0) + 12c(0) 2 )/2T A , 
t = T/T , T A = Aq%/2, and T = Aq 2 B +e /2ttC, where q B 
is the wave vector at the Brillouin zone boundary. 

Near the critical end point y = t /2Ta ~ 0, eq. (B-5) 
can be expanded in terms of y, and we obtain JC(0) — 
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Kq + Kir . The z integral in eq. (B-6) can be carried out 
analytically, and we obtain 



JC(T) 



—7^ f dxx 
Ta Jo 



d+e-i 



logu-^-V(u) .(I3-T) 



where ip is the di-gamma function and u = x e (y + x 2 )/t. 
Rescaling x by t 1 /( 2 + 6 ') j eq. (B-7) becomes 



K,{T). 



/ ds S d+e - 
Jo 



log u — 



1 

2^ 



Here, u = a(y, t, 9)s e + s 2+e with a(y, t, 9) = y/t 2 +<> . As 
discussed by Misawa et a/., 29 eq. (B-8) is proportional to 
if a(y,t,9) -» for t ->• 0. This is understood by 
noting that u — > s 2+e for t — »■ and the integral has no 
t dependence for a small t. 
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